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FINITE ELEMENT AIRCRAFT SIMULATION OF TURBULENCE 


R. E. McFarland 
Ames Research Center 


SUMMARY 


A turbulence model has been developed for realtime aircraft simulation that accommodates 
stochastic turbulence and distributed discrete gusts as a function of the terrain. This model is 
applicable to conventional aircraft, V/STOL aircraft, and disc rotor model helicopter simulations. 
Vehicle angular activity in response to turbulence is computed from geometrical and temporal 
relationships rather than by using the conventional continuum approximations that assume uniform 
gust immersion and low frequency responses. By using techniques similar to those recently 
developed for blade-element rotor models, the angular-rate filters of conventional turbulence models 
are not required. 

The model produces rotational rates as well as air mass translational velocities in response to both 
stochastic and deterministic disturbances, where the discrete gusts and turbulence magnitudes may 
be correlated with significant terrain features or ship models. 

Assuming isotropy, a two-dimensional vertical turbulence field is created. A novel Gaussian 
interpolation technique is used to distribute vertical turbulence on the wing span or lateral rotor disc, 
and this distribution is used to compute roll responses. Air mass velocities are applied at significant 
centers of pressure in the computation of the aircraft’s pitch and roll responses. 


DISCUSSION 


The Simulation of Rotor Blade Element Turbulence (SORBET) model was developed (ref. 1) in 
order to accommodate stochastic and discrete gust velocity components at elements along the rotor 
blades of a helicopter simulation. This model received improved pilot comments over an 
implementation of the conventional MIL SPEC body-fixed formulation (ref. 2) for low-speed and 
low-altitude flight. The introduction of air mass velocities at the blade elements obviates the angular- 
rate filters, because the moment arms to the aerodynamic points of application are implemented. 

This feature has the added benefit that gusts may be input as a function of terrain, so that angular as 
well as translational disturbances may then be correlated with the terrain profile. This is quite 
important for the simulation of low speed flight of various vehicles in the vicinity of buildings, 
ships, etc. 

Recent research has shown that air wake and turbulence modeling require improvement. “A key 
issue with these turbulence models is the way in which they are applied to the aircraft. In all cases 


this was through the aircraft center of gravity. More complex models would require a more effective 
interface with the aircraft” (ref. 3). Insofar as communicating terrain-related information to the 
aircraft model, two options are available. “When airflow gradients are present that are significant 
over distances that are similar to the main dimensions of the aircraft, then it becomes necessary to 
consider introducing forces and moments due to gradients either by ( 1 ) defining gradients at the 
single point and using appropriate derivatives of forces and moments, or (2) calculating local forces 
and moments on different parts of the aircraft and integrating these to obtain total forces” (ref. 4). 

In this material, atmospheric disturbances are separated into stochastic turbulence and discrete gusts. 
Stochastic turbulence is treated in some detail, and differences from the conventional formulation 
are noted in angular responses. Turbulence is developed from the conventional atmospheric 
relationships (ref. 2), except that rotational terms are developed by a distribution of translational 
velocities over the specific vehicle’s effective aerodynamic surfaces using transport delays; a novel 
Gaussian distribution algorithm is used in the development of roll rate as a function of the two- 
dimensional vertical velocity field. A comparison is made to the conventional system. 

Discrete gusts and turbulence magnitudes may be defined from a computational fluid dynamics 
analysis of flow over a terrain profile. It is assumed here that such a database is available; the 
remaining problem is to apply this database to an aircraft simulation model in real time, thereby 
correlating aircraft responses to turbulence as a function of the terrain. 

For a distributed aerodynamic model the appropriate rotational behavior of the vehicle arises 
naturally from force differentials of the various aircraft components, provided only that the 
translational atmospheric disturbances are properly distributed to the appropriate centers of pressure 
(CP) of these components. Such was the case in the SORBET model, which used a rotating blade- 
element model of the UH60 helicopter, with an equal-annuli distribution of aerodynamic centers 
along the blades. However, for total force and moment models, typically using aerodynamic 
coefficients, the vehicle’s rotational behavior in response to discrete gusts requires geometrical 
considerations, which are reviewed herein. 

The realtime transfer of database information to the simulation model, generally including vector 
fields of translational discrete gust velocities and turbulence dispersions, may also include 
atmospheric gradient data. This implies a single-point transfer (probably the vehicle center of 
gravity) for each compute cycle. This formulation appears to be computationally efficient. However, 
atmospheric gradients derived from a grid of translational velocity differences are independent of a 
specific vehicle’s geometry and aerodynamics. In the development of vehicle rotational rates this 
must be considered (see the discussion, ref. 5). Furthermore, using gradients, the implied 
interpolation between data points is likely to smooth significant disturbance parameters unless a very 
large (dense) database is available. Spacing as small as one or two meters have been suggested in 
regions of [flow field] instability (ref. 4). Consider that in a helicopter simulation of landing on a 
ship, “changes in mean airflow in the model appeared too gradual, particularly around the edges of 
the hanger” (ref. 4). “It should be mentioned here that the aerodynamic effects of gust gradients are 
not quite the same as the effects of the airplane angular velocities” (ref. 5). 

In this paper a similar technique to that developed in reference 1 is used to distribute turbulence 
components to the centers of pressure of a conventional aircraft. A model is created called Finite 
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Element Aircraft Simulation of Turbulence (FEAST). With minor modification these techniques 
may be applied to helicopter simulations that use a rotor disc model. 

The compelling reason to use a model like this is to accommodate a component level model of the 
airframe aerodynamics. The distribution of translational turbulence components as well as discrete 
gusts to the various aerodynamic centers would then of itself induce the correct rotational behavior 
of the vehicle. For a coefficient model that does not distribute the aerodynamics, however, we may 
approximate the rotational behavior of airframes by considering only the vehicle-specific moment 
arms to the effective aerodynamic locations. Generating rotational terms for arbitrary aircraft is then 
reduced to a geometry problem. 

The correct moment arms to use in such a formulation would seem to be those that replicate the 
outputs of the conventional turbulence model (ref. 2) in standard flight regimes. It will be seen, 
however, that forcing such a correspondence does not produce realistic moment arms for use in the 
FEAST model, which is developed from elementary mathematical relationships. Rather, when the 
correct separations of centers of pressure are utilized in the FEAST model, the rotational treatment 
in MIL SPEC is brought into question, in particular for the higher frequency regions. This is 
examined in detail. 

In the FEAST model it is maintained that “this idealization of the turbulence field permits 
idealization of the airframe, for calculation of the (generalized) aerodynamic forces and moments, as 
a two-dimensional lifting surface (in the x-y plane) for the vertical component w g and as one- 
dimensional force distributions (along the x axis) for the horizontal components u g and v g " (ref. 5). 
An exponential correlation form is assumed between vertical turbulence velocities. As is shown, this 
form produces considerably more roll activity at low altitudes than does the conventional model. The 
exponential correlation form is certainly not optimal, but improvements in this form are beyond the 
scope of this paper. Referring to scale lengths of a couple of hundred feet, Etkin stated that at “such 
small scales, the variation in gust velocity over the airplane becomes important, and analytical 
methods of some refinement and complexity are indicated” (ref. 6). 

The techniques developed in this paper for stochastic turbulence are applicable to either a gradient or 
distributed finite element model. Terrain related gusts are treated using the specific vehicle 
geometry, and such a formulation should collapse to the gradient method in the mathematical limit, 
providing that the specific vehicle geometry is considered. 

Geometry and Transport Delays 

The FEAST model uses the MIL SPEC translational filters for stochastic turbulence, but does not 
use the angular filters (ref. 2). These filters are replaced by using the geometrical separation of 
centers of pressure, and a novel spanwise correlation factor, as developed from a two-dimensional 
vertical turbulence field. In obviating the pitch and yaw angular filters, the model uses transport 
delays for the emersion rate into the turbulent field. 

The use of a two-dimensional vertical turbulence field permits the development of roll activity as a 
function of the spanwise distribution of velocities, “...when structural modes are significant, the 
exact distribution of the turbulence velocity components over the airframe should be considered ... 
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For purposes of MIL-F-8785B, it is acceptable to consider u g and v g as being one-dimensional 
functions of x, but w g must be considered two-dimensional, a function of both x and y, for 
evaluation of the aerodynamic forces and moments ... The span wise (y) variation of often 
produces significant lateral responses of the airplane and is, therefore, always to be considered, 
except possibly for very slender configurations when approved by the procuring activity.” (Ref. 5.) 

Similar to the techniques used in SORBET, a turbulence “onset line” is defined perpendicular to the 
aircraft body axis in the lateral plane, along which the outputs of the translational turbulence filters 
are computed. For convenience, this line is placed at the fuselage center of pressure. Stochastic 
turbulence components are developed in the body frame, and the onset line is perpendicular to the 
longitudinal body axis. This line is useful in considering the propagation of stochastic turbulence to 
the horizontal and vertical tail centers of pressure as a function of the vehicle velocity. Discrete 
terrain related gusts are later shown to be treated differently from stochastic contributions because 
they are instantaneously applicable at distributed locations. 

Turbulence velocities u M and v M are developed (M=middle line) at the fuselage CP. The vertical 
velocity at the fuselage CP, w M , is produced from two uncorrelated Gaussian inputs given by X and 
Y, as will be shown. Additionally, from these uncorrelated inputs, correlated w r and w ( at the right 
and left wing centers of pressure (along the span) are developed. Only four independent noise 
sources are required, as in the conventional MIL SPEC (ref. 2) formulation, and a distribution occurs 
to the points of application to account for both the geometry and transport lags. 

As the vehicle travels through a frozen field of turbulence, the history of the turbulence encountered 
at a downstream point is equivalent to the history at the fuselage CP, with a time of application 
equivalent to a transport delay. This delay is inversely proportional to the horizontal aerodynamic 
velocity 1 of the aircraft V, and is given by t q = d q /V for the vertical field that causes pitch rate, and 
by t r - d r /V for the horizontal field that causes yaw rate. For example, the lateral turbulence 
velocity applied to the vertical tail’s center of pressure is given by This is merely a 

restatement of Taylor’s frozen field hypothesis, which states that temporal and spatial relationships 
are proportional. 

The yaw rate (due to the tail’s side velocity) is developed by the difference v A/(b-) - v w(0) divided 
by the fuselage-to- vertical-tail separation d r . The sign convention results from the fact that the 
translational turbulence and gust velocities are air mass velocities. The pertinent continuum 
relationships are developed in an early MIL SPEC document (ref. 5) which produced the partial 
derivatives: 


Pg = 



dg ~ 



' “V” is actually the forward velocity of the aircraft, with a lower limit discussed in appendix 1. 
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If the velocity of the air mass is positive downward on the right wing (and absent at the left wing), 
then the rotationally invariant vehicle appears to roll left. If the air mass velocity at the C.G. is 
positive downward (and it is absent at the horizontal tail) then the rotationally invariant vehicle 
appears to be pitching up. If the velocity of the air mass at the C.G. is positive rightward (and it is 
absent at the vertical tail) then the rotationally invariant vehicle appears to be yawing left. The 
translational turbulence and gust velocities of the air mass are subtracted from the vehicle velocities 
in simulation, whereas these rotational terms are added to certain aerodynamic terms (ref. 5). 

The rotational terms may be applied to a C.G.-based aerodynamic model. For a distributed 
aerodynamic model the translational velocities at various aerodynamic centers will induce the proper 
rotations by the development of force differentials. 

For the computation of roll rate, two correlated vertical velocities are produced along the right and 
left wings, at their centers of pressure. The distance d p is the lateral separation of these two points. 
The correlated vertical velocities at the wing CPs are developed from two uncorrelated random 
sequences. As will be shown, these uncorrelated vertical velocities are also used for the computation 
of the fuselage CP vertical velocity w M . 

The generalization of this geometry to a helicopter model requires different points of application. 
These points could include four locations about the hub, typically at radii of about seventy-five 
percent of the rotor radius (ref. 7). In this fashion the rotational gradients due to turbulence on the 
rotor disc may be computed, as well as the influence on the C.G. and tail rotor. 

The Conventional Stochastic Model 

The conventional (MIL SPEC) turbulence model is presented in block diagram form in figure 1. 



Figure 1 . Conventional model. 
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The sums and differences of random variables in this figure are used for convenience. The basis for 
these operations is given in the next section. 

The Dryden form for the longitudinal, lateral, and vertical turbulence filters (that produce u c , v c , w c ) 
are given in Laplace form as follows (refs. 1 and 2): 

fls) _ 

m ' s+V/Z, 




^I3V/kL,,[s + V/V3L,,] 
[s+VILtf 


h{s) 


J3V/jcL„[s + VjSL^\ 

f + WL ,,.] 2 


and the roll, pitch, and yaw rotational filters (that produce p c , q c , r c ) are given by 


lil(s + nV/4b ) 


m(s) 


Ks/Ab 

7+KV/4b 


n{s) 


-ns/3b 

s + nV/3b 


where the quantity b is the wing span. The conventional model’s rotational filters have a history of 
being brought into question. For example, the quantity b has even been replaced by the total aircraft 
length in the pitch and yaw filters (ref. 8). 

In the conventional model, it should be noticed that roll rate is independent of vertical activity at all 
flight conditions (one-dimensional field). The pitch and yaw rates consist of differentiated vertical 
and lateral velocities, using first-order filters with gains that are the inverse of the aerodynamic 
velocity. The z-transform solutions to the above Laplace filters are given in appendix 2. 

The translational turbulence velocities are the same in both the conventional and FEAST 
formulations, where the point of application is the fuselage center of pressure. The rotational 
formulations are different. In particular, FEAST’s roll rate uses certain properties of random 
variables, as examined in the following section. Two correlated vertical velocity sequences are 
developed at the wing centers of pressure, and an exponential form is assumed for the correlation 
between the velocities. The moment arm between these aerodynamic centers is then used to develop 
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roll rate. Pitch and yaw rates in the FEAST model are developed from transport delays acting on the 
vertical and lateral velocity histories. 

Figure 2 summarizes all of the relationships in z-transform notation for the conventional model 
given above, where its outputs have the subscript c. In addition, this figure introduces the FEAST 
model, where its outputs have the subscript s. 

Appendix 1 is provided to show an implementation of the transport delays that are indicated in the 
FEAST paths of this figure. 



Figure 2. The conventional model and the FEAST model. 


Correlating Variables 

In order to create two correlated sequences, two independent, zero mean Gaussian sequences {X} 
and {T} are first developed, both with the same standard deviation, o x = <Jy = <? ■ These may be 
created in realtime simulation by using software 2 such as XNORM©. Two correlated sequences {x} 
and { v} with the same standard deviations o x = a y = <7 (and all other moments) may then be 
created as follows. The correlation coefficient for the sequences {x} and {y} is defined by 


P 


'xy 


(J x O v 


(0<P<1) 


2 The copyrighted Ames Research Center standard simulation software is available to industry. 
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where only direct relationships are permitted. This is accomplished by excluding inverse 
relationships in the selected correlation coefficient, and the lower limit of zero is proper rather than 
negative unity (ref. 9). 

Given the correlation coefficient p as shown above, two symmetrically correlated sequences {x, y} 

may be computed from two uncorrelated sequences {X, Y) with identical statistical properties by the 
functions: 


x ~ 2 (a/ 1 + P + + (V 1 + P ~ 

- v = ^[(Vi+p - + (VT+p + vi r p)i'] 

and their correlation coefficient is p. 

The proof is as follows: E(XY ) = 0 for independent variables, and E^X 2 ) = e[y 2 ) = a 2 . The 
expected values for the dependent variables are given by: 

£ {* 2 } = \ E {[(^ + + (^p - ^ z p)yf] 

= ^g{ (v/f+P + a/ 1_ p) X 2 + (v/ 1+ P ~ y/ 1_ p) 2}/2 } 

= ^ 2 [( 1 + P)+(1-P)+(1 + P)+(1-P)] = o- 2 

E {y 2 } = | £ {[(Vi+p - + (Vi+p + Vkp)^] 2 } 

- 4 £ |(V 1 + P ~ V 1- p) + (y/ 1 + p + v/ 1_ p) ^ 2 | 

= | cr2 [( 1 +p)+( 1 -p)+0+p)+(i-p)] = tf 2 

£{yy} = ^{[(V^ + V^J^+fV^-V^J^V^-V^^+fVi+P+Vi-p) 1 "]} 

= {g{(V^ + V^)* 2 + (V^- + V^)*' 2 } 

= | fJ ' 2 [(Vl + P + v/i-p)(vA+P - v/i-p)] = <* 2 P 

Thus 


cr 2 C7,cr v 


(0 < p < 1) 
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In the FEAST model this result is used to produce correlated turbulence velocities at the wing’s 
spanwise centers of pressure by relating 3 the vertical velocities x => w r and y => w? as linearly 
filtered random variables. These are computed from the uncorrelated X => w R and Y => w L as 
shown below. An exponential correlation form is used: 


£[w(r 0 )w(ro + r)] = a 2 e & 

where for this model L is the characteristic length 4 of the vertical component of turbulence and r is 
the distance between the two points of interest (the difference between the wing centers of pressure, 
d p ). Thus, the correlation coefficient becomes 


P 



and the criterion for only positive correlation is satisfied due to the lower limit of the MIL SPEC 
definition (ref. 2) for low-altitude turbulence models 


10 

h <\0 ft 

h 

10 ft<h< 1000 ft 

1000 

h> 1000 ft 


where h is the vehicle altitude. For a center of pressure difference d p = 32.17 ft, the correlation co- 
efficient would be as shown in figure 3. 



Figure 3. Wing CP correlation coefficient. 


' The random sequences first pass through the H(z) atmospheric filter. Gaussian variables that pass through linear filters 
remain Gaussian variables. 

4 The characteristic length is the spatial separation of two points where the correlation of velocity components is 
assumed to vanish. 
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The vertical turbulence velocity at the fuselage CP may also be computed from the uncorrelated 
variables. This quantity vv M results from the correlation algorithm by setting p = 1 , a condition 
called for when the distance between the two points of interest vanishes (r = 0). Using the uncorre- 
lated quantities X and Y we see that this produces the linear combination for the centerline vertical 
turbulence velocity 



W M = Jjl w L + Wfl] 


which is the Gaussian centroid, and all statistical moments are retained, as in the case of the 
turbulence at the wing span positions. X and Y are filtered to produce w R and w L . 

To complete the picture for characteristic lengths, the MIL SPEC restricted lower limits for the 
horizontal components are also used, and it is noted (ref. 2) that the off-axis dispersions are deter- 
mined once the vertical component a w is selected. 


Ml 


h 

(0.177 + 0.000823/j) 1 ' 2 


(10 ft < h < 1000 ft) 


a 


u 


C7, 




(0.177 + 0.000823/?) 0 ' 4 


(10 ft < h < 1000 ft) 


Using the relationships of this section, the spanwise center-of-pressure vertical velocities uy(z) and 
Wf(z) for the right and left wings, and the midline vertical velocity w^(z), are created from two 
independent noise sources. An example time history of these velocities is presented in figure 4, 
where the correlation coefficient used was p — 0.6525. 



Figure 4. Vertical turbulence at wing and fuselage CPs. 


The velocities in the above figure were created at an altitude of 40 feet, an aerodynamic velocity 
of 5 ft/sec, and d p = 17.08 ft. The CP vertical velocity histories were used to produce the roll rate 
given in figure 5. 
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Figure 5. Resultant roll rate history. 

Note that correlation between the vertical velocities at the wing CP locations is only a statistical 
measure. It does not provide any indication of the direction of roll with respect to the fuselage 
vertical velocity. Indeed, roll is completely uncorrelated with the fuselage vertical velocity, as the 
following expected value shows: 

E{p s w M } « E{[w l - w R \w L + w/f]} = a 2 -a 2 = 0 


Vertical Velocity 

Consider the vertical velocities W L and w R , which are uncorrelated Gaussian random variables with 
the same dispersion <7, and zero mean value. Such variables with identical statistics have the 
following properties: 


e{w l w r } = 0 


£{[W l + W*] 2 ] 

■ = E{wl + w*\ = 2 a 2 

e{[w l -w r ] 2 } 

■ = e{ w l + w r\ = 


E{[w L + w R \w L -w R ]} = E{w 2 -w 2 } = 0 


The expected values of both the sums and differences produce new random variables that have twice 
the variance of the originals. Hence, we may create new random variables with the same dispersion 
a by dividing by the square root of two. That is, in terms of discrete, unity-variance random 
variables 
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and from the above properties the sum represented by ri 5 (z), the Gaussian centroid, is uncorrelated 

with the difference represented by t] 6 (z). For practical purposes the random turbulence vertical 
velocity may be created in z-transform notation from the sum 


w(z) 


[K h( 7 \ 

'm(z)+ih(z)~ 

y j 

V2 J 


a w ^ff(z)rj 5 (z) 


where the system function for vertical velocity in Laplace notation is 


vv(i-) = O w h(s) 


I 3V (s + V/V 3L^j 

(.v + V/L w ) 2 


The autospectrum or power spectral density (one-sided PSD) of the Dryden form of the vertical 
response is given by (ref. 2) 


3 Vcr 2 (ft) 2 + P 2 /3 l 2 ,) 
kL w (co 2 + V 2 // 2 ,) 2 


“The total area under the autospectral density function is the variance of the data plus the square of 
the mean value of the data” (ref. 9). Since the mean value of the random input is zero, the variance 
of this autospectrum is given by the frequency integral from zero to infinity. 


ot = 


1 


<£> w (co)dco = 


3V<t- 

ttLvj. 


K 


3VC7 2 

kU, 


r 




CO 2 + 


y2 l A 2 v (co 2 + V 2 /lZ,) 2 


2(V/A v .) 


- |(v/Av) 2 


K 


4 (V/L^) 3 


\dco 


= ot 


Rotational Variances 

In appendix 3, the conventional formulation’s rotational relationships are reviewed and the rotational 
variances are shown to be 



0.4<7 2 ;r 
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7^ = 

’ qc 


<F 2 w 2f3«t 2L + 2 

a w ,7r ^ 4& 


32b 2 (^ + l 
V 4b 


a ^i^t +2 ) 


Also, in appendix 4 the FEAST model’s rotational relationships are developed. The rotational 
variances are shown to be 


2al(l-p) 
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Equating Variances 

Equating the variances between the conventional and FEAST models would seem to be an excellent 
technique for quality control. In pitch and yaw the sharper responses of the FEAST model could be 
used to justify the lower-frequency continuum responses of the conventional model, and in roll the 
influence of the vertical correlation relationships could be examined to determine the altitude where 
behavior was similar. 


Unfortunately, as determined in this section, the technique of equating the variances produces 
completely unrealistic CP separation distances. 


For comparison of the two formulations with equal variances, the distance between the wing centers 
of pressure for use in the FEAST model may be determined if the span ( b) and the characteristic 
length ( L*,) are known. Solving for d p in the roll dimension, this requires the solution to the 
following transcendental equation, obtained by equating o^ c — 


0.2 n\ 


4 b ) 


- d p/ L 
1 -e /Xj>< ’ 

(dp/L*,} 
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Similarly, equating cr qc = a qs and <jy C = also produces transcendental equations 

3 + 2{4b/KL^) _ 2 + [(^/L»,)-2]c 

2(4/?/ nL, v )[ 1 + (4b/ nLn )] 2 K/'Av) 2 

3 + 2(3£/ttL 1 ,) _ 2 + [(t/ r /L,.)-2]c-^ /L ’ 

2{3b/7TL,.)[\ + {3b/KL v )f {d r /Q 2 

which may be used to solve for the other two CP differences. These three equations are independent 
of both variances and vehicle velocity. They are easily solved using Newton’s iteration (ref. 10). 
Very few cycles are required. 

For i»=32. 17 ft and 7^=250 ft (7^=791.48 ft) Newton’s iteration produces the values given in 
figure 6 for the center-of-pressure separation distances d p , d q , and d r . 



Figure 6. Computing the required CP differences using equal variances. 

As shown in figure 7, the separation distances vary somewhat when equal variances are used 
(£=32.17 ft). 
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o 200 400 600 800 1000 
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Figure 7. Separation distances vs. altitude for equal variances. 

The CP separation distances for pitch and yaw are fairly constant for high-altitude flight, although 
for spectral correspondence with the conventional formulation the required CP values for the FEAST 
model are too high (not physically credible) for typical aircraft geometry. The required separation 
distance for roll correspondence is also too high, except that it appears to approach a reasonable 
value for high-altitude flight. 

Although the variances are independent of vehicle velocity, the autospectra are not. Selecting 
both a low and moderate velocity case for illustration, the following six figures show the spectral 
differences between the conventional and FEAST models, where the CP separation distances 
have been computed from the equal-variance relationships. In these examples where L, v = 250 ft 
= 791.48 ft), the selected value for o w = 1 ft/sec, which produces a v = 1.468 ft/sec. The 
first three figures are for a velocity of 5 ft/sec. 



Frequency - rad/sec 


Figure 8. Roll autospectra (V = 5 ft/sec) using equal variances. 
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b = 32.17 ft, d q = 81.07 ft 
L = 250 ft, o =.0110 rad/sec 

w q 


Frequency - rad/sec 

Figure 9. Pitch autospectra (V = 5 ft/sec) using equal variances. 


^rsM 



Frequency - rad/sec 

Figure 10. Yaw autospectra (V = 5 ft/sec) using equal variances. 
The second set of three figures is for a velocity of 100 ft/sec. 


b = 32.17 ft, d p = 33.39 ft 
L = 250 ft, a = .01 50 rad/sec 

w p 


Frequency - rad/sec 

Figure 11. Roll autospectra (V = 100 ft/sec) using equal variances 








Figure 12. Pitch autospectra (V = 100 ft/sec) using equal variances. 



Using Reasonable CP Separation Distances 

For an aircraft with a span of 32.17 ft, the computed d p , d q , and d r required to obtain equal 
variances are clearly too large. An example of an aircraft that has this particular wing span is the 
XV- 15. Data is available for the center-of-pressure distances for this aircraft. Considering the CP of 
the fuselage given by 24.416 ft, for the XV- 15 the differences in these distances are 

d p = 17.08 ft 

d q = 46.666-24.416 = 22.25 ft 

d r = 47.501-24.416 = 23.085 ft 

b = 32.17 ft 

We may investigate the changes in FEAST’s autospectra using these values. As above, both low and 
moderate speed cases are shown in the following six figures. The first three figures are for a velocity 
of 5 ft/sec. 
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Yaw Autospectra Pitch Autospectra Roll Autospectra 



Frequency - rad/sec 

Figure 14. Roll autospectra (V = 5 ft/sec) using CP distances. 



Frequency - rad/sec 

Figure 15. Pitch autospectra (V = 5 ft/sec) using CP distances 




^ ^rcM j 

b = 32.17 ft, d r = 23.09 ft, = 250 ft 
a =.0112 rad/sec, a =.0186 rad/sec 


10' 3 10' 2 10' 1 10° 

Frequency - rad/sec 

Figure 16. Yaw autospectra (V = 5 ft/sec) using CP distances. 



Yaw Autospectra Pitch Autospectra Roll Autospectra 


The second set of three figures is for a velocity of 100 ft/sec. 



Figure 17. Roll autospectra (V = 100 ft/sec) using CP distances. 



Figure 18. Pitch autospectra (V = 100 ft/sec) using CP distances. 
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From these figures we see that the pitch and yaw spectra in the FEAST formulation are extended in 
the frequency domain, and they have larger dispersions because of this. The roll spectrum has less 
bandwidth, but a larger dispersion. These phenomena are well within a pilot’s bandwidth, and would 
thus impact simulation results. The resultant rotational dispersions are compared to the conventional 
dispersions in the following three figures. a c is the conventional model’s dispersion and a s is the 
FEAST model’s dispersion. 



Figure 20. Roll dispersion. 



Figure 21. Pitch dispersion. 



Figure 22 . Yaw dispersion. 
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Step Responses through the Turbulence Models 

The conventional rotational turbulence model uses only the wing span parameter, where the value 
b = 32.17 ft corresponds to the XV- 15 aircraft. When the FEAST model is used, and it is 
statistically required to produce the same angular dispersions as the conventional model, un- 
realistic CP differences were computed, given by d p = 33.39 ft, d q =81.07 ft, and d r = 61.39 ft. 
For the XV- 15 aircraft, the smaller values of d p = 17.08 ft, d q = 22.25 ft, and d r = 23.09 ft are 
known, where the span length is also given by b = 32.17 ft. These three configurations produce 
three different step responses, given in the following three figures. These responses were created 
at an altitude of 250 feet, and a vehicle velocity of 5 ft/sec. 



Figure 23. Roll rate, step input. 



Figure 24. Pitch rate, step input. 
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Figure 25. Yaw rate, step input. 


From an examination of these figures the conventional model is shown to be a low frequency 
approximation to the behavior computed in the FEAST model. The artifact of setting the dispersions 
equal and then computing the required CP separation distances produces somewhat similar pitch and 
yaw responses, but FEAST’s roll response is larger. 

In all cases it is seen that the conventional model has less rotational activity than the FEAST model, 
and this is especially true if the correct CP distances are used. 

Although step responses are not particularly revealing in terms of responses to random data, these 
figures do indicate that the conventional model’s rotational responses may only be valid for low 
frequency behavior. 


Time Response Comparisons 

In the following material the XV- 15 aircraft’s CP separation distances are used in the FEAST model, 
and the XV-15 aircraft’s span length is used for the conventional model. The velocity is 5 ft/sec, and 
time histories are presented for two different altitudes. Identical random sequences were used to 
drive both formulations. 
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It should be noticed that the scales used in these figures for the low altitude cases are much larger 
than those used for the high altitude cases. 

The FEAST model roll behavior tends to approach the conventional model roll behavior for high 
altitude flight. Pitch and yaw responses, however, are generally larger for the FEAST model. 

In reference 1, a technique was used to distribute random turbulence over the UH60 helicopter’s 
rotor disc. It ensured that the dispersion was identical at all rotor blade elements by using an 
algorithm that is equivalent to the one developed herein for roll activity. It did not assume any 
particular correlation form. Because the UH60 simulation is a distributed aerodynamic model, 
individual element forces were computed from extensive aerodynamic tables, and moments were 
induced from the rotating frame geometrical relationships. The results of this SORBET simulation 
model, which had extensive piloted testing, may be used to reveal a feature in common with the 
FEAST model. Although pilot opinion was generally favorable for the SORBET model, the 
following quote indicates that the turbulence dispersion should be reduced for low velocity flight 
(which implies low altitude). From reference 1: 

“Pilots generally agreed that the turbulence magnitude should be a function of velocity. 

Above 40 knots, the subjective values for light and moderate turbulence were selected as 
(7 W = 5 ft/sec and 8 ft/sec, respectively. Near hover, pilots selected standard deviations that 
were half these values. The reason for this is unknown, and should be investigated further.” 

The FEAST model’s increased angular dispersions of figures 20-22 are thus an approximation to 
what would occur if a distributed aerodynamic model were used, and the objectionable magnitudes 
may be isolated to the translational spectra at low velocities and altitudes. The integrals of these 
spectra (the translational variances), however, are selected, and they are independent of velocity and 
altitude. Hence, it is suggested that the translational variances themselves should have some 
functional decrease with altitude (in the low altitude region). This is strongly suggested (ref. 2) in a 
graph of <r w versus altitude, which is unfortunately applicable to “medium/high altitude only.” 
Clearly, further research is required to determine the proper functionality for low altitude flight. 

Discrete Gust Contributions 

The literature is replete with papers dealing with mapping the environment around ships (refs. 1 1-13). 
Also, current projects at Ames Research Center include terrain mapping of significant features of a 
visual data base to correlate gusts to these features using computational fluid dynamics (ref. 14). The 
mapping functions are exceedingly complex; these complexities are avoided here by an example using 
a simple terrain related gust model. The FEAST model produces complex atmospheric/vehicle motion 
in response to a simple terrain related model. 

Gusts are assumed to be created in a frame of reference corresponding to the visual model, and 
transformed to the instantaneous distributed centers of pressure. Terrain correlated gusts may be 
formulated to produce representative angular activity if they are applied at the pertinent CPs. A 
system that uses this geometry requires the transmittal of seven gust components from the 
atmospheric model to the aircraft model. These consist of one u-gust, two v-gusts, and four w-gusts. 
The fuselage CP and the CPs of the tail and wing are the pertinent points. The quantities required 
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each cycle time are: u g ( 0), v g (0), v g (d r ), Wg(0), w g {+\dp), w g {-\dp} and w g {d q ^. The center- 
of-gravity translational discrete gusts are then the values with the arguments of zero (the CG is very 
close to the CP). The rotational contributions are then given by the equations 

_ H g(~l^) ~ H g( + 2 dp) 


w g {0)-w g (d q ) 
q g d q 

v g{d r )-v g { 0) 

“ d r 

An example using these discrete gust equations is given in appendix 5. 

The fuselage CP vertical discrete gust could be approximated by the average of the wing CP gusts, 
which would reduce the required number of gust transmitted quantities to six. 

w g{~^p) + w g{ + 2^p) 
w g ~ 2 

The total contributions to the atmospheric model are given by sums of the stochastic turbulence and 
discrete gust velocities, such as u a (z) = u s (z) + u g (z) and p a (z) = p s (z) + P g (z)- 


CONCLUSIONS 


A generic turbulence and gust model has been developed for aircraft simulation, where both the 
translational and rotational components are computed. This model, called FEAST, accommodates 
stochastic turbulence and distributed discrete gusts as a function of the terrain, while the required 
aircraft-specific parameters are related to the vehicle geometry. FEAST delivers rotational stochastic 
components that are extended in the frequency domain beyond those delivered by the conventional 
model. 

New discrete techniques for the computation of rotational rates are introduced. A correlation 
coefficient based upon the physical separation of aerodynamic centers, and a two-dimensional 
vertical turbulence field are used in the computation of roll rate. Filtered random noise sources are 
propagated in time by use of interpolated transport delays. For the computation of pitch rate, the 
Gaussian centroid vertical velocity of the field influences the centers of pressure of the fuselage and 
the horizontal. tail as a function of vehicle velocity. Similarly, the vertical tail’s center of pressure is 
used for the computation of yaw rate. 
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The FEAST model should be useful for low speed and low altitude flight simulation. It uses the 
specific vehicle geometry in the computation of rotational responses to turbulence, and gusts as a 
function of the terrain. The state-space techniques given for the discrete implementation of this 
model are computationally efficient, and deliver stable outputs for all stable inputs regardless of the 
required cycle time. 
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APPENDIX 1 


Interpolation and the Circular Buffers 

The pitch center-of-pressure difference d q may be larger than the yaw center of pressure difference 
d r , as is assumed in this appendix. Regardless, the larger value is used to compute the minimum 
aerodynamic velocity that must be imposed in order to retain all of the required vertical and lateral 
velocity data in buffers of size N. In order to accomplish this we have 

d a 

V > 

V {N-\)T 


The velocity used throughout this paper is assumed to conform to this inequality, resulting in a 
minimum velocity of a couple of feet per second using reasonable buffer length and cycle time 
parameters. We may then develop an integer index k q with the following limits 5 : 


0 < k q 



< N - 1 


with remainder P q given by 


0 < P q 



k q < 1 


The lower velocity limit assures us that if the first ( n = 0) cell of a buffer contained the current time 
point’s fuselage CP data, then we could interpolate within the buffer to find the horizontal tail’s CP 
vertical velocity data by use of 


H’ 


VT 


[l - P q ]w(k q ) + P q w(k q + l) 


and find the required end-point values, regardless of the aerodynamic velocity (above the minimum). 
Similarly, if the yaw CP difference is smaller than the pitch CP difference, then the same minimum 
velocity is used such that the lateral inequalities become 


0 < k r = 

0 < P r = 



N- 1 


< 1 


and the vertical tail’s CP lateral velocity data may be obtained from 


5 Using the “floor” or “ intier ” operation (least integer). 
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\vf) = t 1 “ PrK k r) + PA k r + 1) 


The buffers are initially full of velocity data. This may be accomplished by a negative-indexing 
initializing scheme. Using the filtered random inputs, first the last cell ( N - 1) is filled followed by 
(N — 2), etc., until the index becomes (0). In simulation this requires N passes through the 
turbulence equations (initial condition mode) before realtime operations are begun (operate mode). 
The first storage cell in realtime is then also 6 N- 1 , and thereafter, a decrement continues each cycle 
time for the entire run, with a modulus of N. 

n => (n-l)Mod(N) = (n-l)-N -- - - - 


The “right arrow” is here read “becomes”; at the transition point this produces 


0 


(0 - l)Mod(N) 


(O-l)-Al 


0-1 

N 


1 -N\ 



N- 1 


effectively completing the circle. The circular buffers each require only one store operation every 
cycle time, into the cell designated as n . The next higher indexed cell (ModN) is always the cell 
from the previous deposit, so this scheme is sufficient for linear interpolation of the delayed vertical 
and lateral velocity values. 

From the above preliminaries, the entire procedure is given below in twelve steps for depositing the 
current filter outputs into the circular buffer (where these values are the fuselage CP velocities), and 
retrieving the vertical and lateral velocity at the horizontal and vertical tail centers of pressure. 


Step (1) 

n => 

(n-l 

\N~ 

n > 1 
1 n < 1 

Step (2) 

Deposit current filter 

Step (3) 

k q = 

1 


Step (4) 

A = 

dg 

VT~ 

- b 

Step (5) 

m q = 

n + k q - N 


n + k q 

~N~ 


6 Not really required as long as the buffer is full and the correct index is retained in transferring to operate mode. 
Initializing just requires at least N cycles. 
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Step (6) 

Step (7) 
Step (8) 
Step (9) 

Step (10) 
Step (1 1) 
Step (12) 


* _ f /n <? + 1 m q <N ~\ 

" z< ? ~ [0 m q -N-\ 

k = — 

VT 



W\ 


(d^A 

k vt) 


J4l 

\VT 


n + k r - N\ 


n + k r 

N 


jm r + 1 m r < N - 1 
(0 m r = N — 1 


= [l~ Pq]v(™ q ) + Pqv(m*) 

= [l - P r ]v(m r )+ P r v(m* ) 



APPENDIX 2 


z-Transform Solutions 

For discrete implementation of the translational filters using either of the two models a zero-order 
hold formulation is typically used, as shown by the three z-transform relationships: 


F(z) = Z 


l-e 


-sT\ 


f(s)\ = 


_ 2L, ( 


nV 


~T t ^ 

l-e L “ 

[Z-e 


G(z) 


fl- e sT\ 





m = z 


f\-e- sT ^ 





where the coefficients are: 


~T 


VT -~ T 


C, = l-e L v + U?>-l)^j-e L ' 


C? — — e 


U 


~T 


+ (V3-l)^ 


c 


l-e 


-y~T 


(V3-1) 


VT -U T 


Lw 


Ca 




l_ f ,‘i. r + (V3-l)tZ 


The difference equations that result from the z-transforms are driven by zero-mean, unity variance 
discrete Gaussian noise sequences r\ k (k=l, 2, 3, 4) with amplitude range {-3.45, 3.45}, which 
preserves 99.94% of the continuous total probability density function. These derivations assume the 
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utilization of the realtime XNORM© Gaussian generation software, with zero mean and unity 
standard deviation. 

The conventional model’s rotational terms are developed below to make a comparison with angular 
outputs of the FEAST model. The rotational terms require three z-transforms. Whereas the roll L(z) 
term uses the zero-order hold formulation, as was used in the development of the translational z- 
transforms, the pitch M(z) and yaw N(z) formulations use the triangular hold, because they are 
computed, in sequence, after functions that use the zero-order hold. This preserves the proper 
input/output phase relationships between all variables. 



For discrete implementation a normalizing factor is required to preserve the total power integrals 
(cr 2 ). Where the 77,(2) are unity-variance, zero-mean Gaussian variables, the translational z- 
transforms are given by 

“c(z) = e u ^ F (z) 7 h(z) 

v c (z) = o v ^G(z)ti 2 (z) 

w c (z) = = cr w ^^H(z)[Ti 4 (z) + T] 3 (z)] 

and, with a sampling switch between translational and rotational filters, the rotational rates become: 
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Pc(z) = O w ,^Y L(z)tj 6 (z ) = o w ^^L(z)[v4(z)-rh(z)} 

q c (z) = (T w ^H(z)M(z)rj 5 (z) = a w ^;H(z)M(z)[ri 4 (z) + ri 3 (z)] = M{z)w c (z) 
r c (z) = ct v ^G(z)N(z)tJ2(z) = N(z)v c (z) 

Note that these rotational filters are not used in the FEAST model. 

Difference equations are easily created from these z-transforms. For example, the longitudinal 
turbulence velocity becomes 


u c (n + l) = e L “ u c {n) + o u 

State space solutions to the turbulence filter equations are preferred because they produce stable 
outputs from all stable inputs, independent of the size of the cycle time parameter, and the inputs to 
these filters have considerable frequency content. Additionally, the zero-order hold solution 
produces an exact solution to a step input, and the triangular hold produces an exact solution to a 
linear input. Note that the discrete coefficients must be recomputed whenever velocity or 
characteristic length change. Due to the fact that these are “slowly varying” quantities, however, a 
practical implementation technique is to distribute the workload of the coefficient-update operation 
over a number of computer cycles. In this way the workload becomes competitive with that of less 
robust algorithms, and the entire coefficient set is updated in a few computer cycles. This technique 
has been used at Ames Research Center for 20 years. 
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APPENDIX 3 


Conventional Model Rotational Relationships 

The system functions, autospectra, and variance are given for the conventional model. 

Conventional Model Roll Relationships 

The conventional (MIL SPEC) formulation’s roll response may be written in terms of the 
independent unity-variance random variable ^(z) 

Pc(z) = <? w ^P(zMz) 


For this conventional formulation the roll system function is 

o w V(X 8 V (k/ 


Pci*) = 


The roll autospectrum given by 




L%{s + nV/4b) 


0.8 alyijtj4bf^ 


L$[(0 2 +(7CV/ 4b) 2 ] 

and its variance is computed by the frequency integral from zero to infinity 

0.8cr^ V(^/4fe)^ f dco 

Jo 


r 2 = 

' pc 


i. 


<s> pc m<0 - £ J 0<j) 2 +(nv/4b y 


O.SalV f n Y /j 
= i% U b) 


K 


2(nV/4b) 


0.4<r;|,7r ( n_ Y* 

1% V^b) 


Conventional Model Pitch and Yaw Relationships 

Where a sampling switch is assumed between the linear and rotational filters, pitch and yaw for 
conventional model are given in z-transform notation by 



q c (z) = w{z)M(z)ri 5 {z) = (T w ^yH(z)M(z)t] 5 (z) 
r c (z) = v(z)N(z)t] 2 (z) = <J v ^G(z)N(z)t} 2 (z) 


The system functions for these rates may be written 


* w = aMs)m(s) = 


r c (s) = cr,v(s)n(s ) = - 


° vK ' 3V (s+V v3L,.).v 


3 b yJtLv 


(.v + 7tV/3b)(s + V/L,) 


The autospectra are given by: 






^^(0 2 ((0 2 + V 2 / 3Z 2 ,) 

16 b^, \ ' 1 

(to 2 + n 2 V 2 !\6b 2 ^{(o 2 + V 2 / L 2 ^ 

°vVk 2 / 2 


3b 2 Ly 


(0 2 {o) 2 + V 2 /3Z^) 


((0 2 +7t 2 V 2 /9b 2 )((0 2 + V 2 /l 2 ,y 


and from these expressions the variances are: 


f 


tfc = I <& qc {a>)d(a 


J, 


<& rc {(0)d(d = 


aW(^ + 2 


v 46 


32b 2 ( KL "’ 


46 

f ?rA. 


+ 1 




186' 


(riv 
{ 3 b 


+ 1 
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APPENDIX 4 


FEAST Model Rotational Relationships 

The system functions, autospectra, and variance are given for the FEAST model. 

FEAST Model Roll Relationships 

Using the properties of random variables, an expression for roll rate in the FEAST model may be 
developed. Where d p is the difference between the centers of pressure for the left and right wings, 
vertical velocities applied at these points produce a roll rate given by 

W( - w r 

Ps - 2 

Up 

These velocities at the centers of pressure are correlated. However, using the previously developed 
correlation coefficient the roll rate may be expressed as a function of the two uncorrelated vertical 

velocities wi and w R 


Ps = 


_ W( — W r 


-p[w L -Wft] 


Hence, FEAST’s roll rate in z-transform notation becomes 


p,(z) = = SzdMElM H(zMz) 


and the system function for roll rate is given by 

PsO ) = 


jW( s+v /^ L ^) 
d D W{S) d p inK (s+v/ r) 2 


FEAST’s roll autospectrum is thus a function of the vertical autospectrum 


<1 vM 


2(1 -P) 

4 




6ojV(j - p)(m 2 + V 2 /3 ij) 
JtL w d 2 p ((0 2 + V 2 lll) 1 


although roll response is uncorrelated with vertical response. The roll variance is similarly computed 
by a frequency integral from zero to infinity. 


ps 


f 




2(1 -P) 


r*. 

Jo 


{(o)do) 


2o 2 (l-p) 
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FEAST Model Pitch and Yaw Relationships 


It can be assumed that w(z) is created at the fuselage center of pressure without loss of generality if 
transport delays are then imposed, which are equivalent to the distances to the horizontal stabilizer 
( d q , for pitch rate) and vertical stabilizer (d r , for yaw rate). The time shifts are a function of the 
vehicle velocity 



which, in discrete simulation, are equivalent to selecting and interpolating filtered velocity values 
from circular buffers as shown in appendix 1. The pitch and yaw rates in z-transform notation are 
then given by. 


<7j(z) 

»i(z) 


w{z) 


) dr 


These produce the system functions 


K_ 

T 

( 


-Vd 


H(z)Vs(z) 


-iRzfeW 


q s (s) 

r s(s) 




a w h(s) 



& V g(s) 


dq{s + V/ Lh ,) 2 

dris+VIL,) 2 


which may be written generically 


C,(s) = 


d(s + a ) 2 


with a generic autospectrum given by 

3ao 2 (2 - e _;<UT - e^ m + a 2 j 3) 




7td 2 [co 2 + a 2 ) 


6a<7 2 (\ - coscur)(ty 2 + a 2 / 3) 

xd 2 (o> 2 +a 2 f 


The individual autospectra thus become 
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<M<y) 


6 ( 7 2 V(l -cosft*^/v)(ft) 2 + V 2 / 3 /£.) 
kL, , d 2 (a ) 2 + V 2 /L 2 v f 


<M") 


6 < 7 2 V(\ - cos ( 0 d r jV){(D 2 + V 2 / 3 Z^) 
nL,,d 2 (( 0 2 + V 2 /! 2 ,) 2 


The generic variance is given by the integral 


& = 


i 


O s (a))da) = 


6 flg 2 r 

Jc 


(1 - cos<wt) 


2 2 
— n 

3 “ 


® +fl (® 2 + a 2 ) 




6 a cr 
7itf 2 


f 


2 2 
——a 


C0 2 + a 2 3 


(ft ) 2 + a 2 ) 


COSft)T 2 2 COSft)T 
+ — a 


2 ft) 2 + a 2 3 


(ft ) 2 + a 2 )‘ 


Wft) 


6 acJ“ j 7 T 2 2 
~nd 2 ~\ 2« 3 fl 


(ar- 2 )e aT +2 
J 2 


K 

ne- ax 2 2 
~ 

-^-(l + aT)f? aT l] 

.4 a 3 . 

2 a 3 

. 4 a 3 V ' JJ 


From this expression the FEAST model’s pitch and yaw variances are given by 
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APPENDIX 5 


Discrete Gust Example 

To illustrate the influence of the distributed points of application in relationship to a terrain model, a 
simplified canyon (or ship deck) is constructed. This model uses two separate terrain altitudes, in a 
Cartesian plane, as shown in figure 27. 

y 


Figure 27. Terrain elevation. 

An aircraft travels at a constant speed (V), at an azimuth angle (y/), and at a constant altitude. The 
altitude of the aircraft above the plateau is 50 ft, and over the canyon it is 100 ft. For illustrative 
purposes, neither the aircraft trajectory nor attitude is influenced by the gusts. The disturbance model 
is constructed on this plane such that a vertical (upward) gust of 5 ft/sec is encountered when the 
aircraft arrives at the canyon wall. This occurs at Xj(t) = 0, where / = ( F , HT, VT, RW, LW ) denoting 
the fuselage CP, the horizontal tail CP, the vertical tail CP, the right wing CP, and the left wing CP. 
As a function of time (t) the points of application are given by 

x F (t) = x F (0 ) + Vtsiny/ 

xht( 0 = x F (t) - d q sin y/ 

xvr(t) = *f(0 - d r sin y/ 

*Rw(t) = x F {t) + {d p cosy/ 

x iw{t) = x f ( 0 - {d p cosy/' 

y j(t) values are not required in this example because the y-gust is assumed to be a function of the w- 
gust, which is only a function of *,(/). For convenience, the vertical gust decays linearly to zero 
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when Xj(t) = 100 ft. Also, a positive y-gust of one-half the magnitude of the vertical gust is assumed 
to exist that is independent of the y ( (f) value. 


fo 



10 


- 100A 
100 ) 


Xj < 0 

0 < Xi < 100 ft 
x, > 100 ft 


l g(i) 


-w g{i) cosy/ 


v g(i) = 2 W ^ Sin¥ 


The discrete gust profile given above leads to asymmetric responses when the canyon is approached 
at an angle i/r = 45°. In order to investigate these responses the initial condition x F (0) = -17.68 ft is 
used, where y = 5 ft/sec. The rate of closure on the canyon becomes 3.536 ft/sec, so that it takes 
exactly 5 seconds for the fuselage CP to arrive at the canyon. Due to the geometry, (1) the right wing 
arrives first at t =3.29 sec, (2) the fuselage CP arrives at t = 5.00 sec, (3) the left wing arrives at 
t=6.71 sec, (4) the horizontal tail arrives at t=9.45 sec, and (5) the vertical tail arrives at t=9.62 sec. 
The vertical gust velocity that occurs at each of these aircraft locations is given as a function of time 
in figure 28. 



Figure 28. Vertical gusts at aircraft locations. 


The gusts influence the fuselage CP in translational space. The onset of these gusts thus occur at 
t = 5 sec, and they are as shown in figure 29, where the u-gust and v-gust are half the magnitude 
of the w-gust. 
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Figure 29. Translational gusts at the center of gravity. 


Induced rotations occur as the various aerodynamic centers become immersed in the gusts. Since 
the right wing arrives at the canyon at t = 3.29 sec, the roll response begins at this time. Note that 
partial recovery occurs at t = 6.71 sec, when the left wing also becomes immersed in the gust field 
The differential gradient occurs due to the different vertical gust magnitudes, out to a point where 
neither wing remains in the gust field. The pitch and yaw responses occur almost simultaneously 
because d q ~ d r . These responses are shown in figure 30. 



Figure 30. Rotational response to discrete gusts. 

From this example we see that the FEAST model produces atmospheric variations based upon the 
specific vehicle geometry in response to terrain-related discrete gusts. This is important in evaluating 
responses of different aircraft to the same terrain related (or ship) turbulence and gust profile. 

Refinements to the FEAST model await the conversion of an area of CGI-terrain data to an 
appropriate form for analysis using computational fluid dynamics, and the application of the 
resultant flow-field data to various aircraft simulations. 
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